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ABSTRACT 

Aims. The goal of this work is to measure the evolution of the Galaxy Stellar Mass Function and of the resulting Stellar Mass 
Density up to redshift ~ 4, in order to study the assembly of massive galaxies in the high redshift Universe. 
Methods. We have used the GOODS-MUSIC catalog, containing ~3000 Ks-selected galaxies with multi-wavelength coverage 
extending from the U band to the Spitzer 8/im band, of which 27% have spectroscopic redshifts and the remaining fraction 
have accurate photometric redshifts. On this sample we have applied a standard fitting procedure to measure stellar masses. 
We compute the Galaxy Stellar Mass Function and the resulting Stellar Mass Density up to redshift ~ 4, taking into proper 
account the biases and incompleteness effects. 

Results. Within the well known trend of global decline of the Stellar Mass Density with redshift, we show that the decline of 
the more massive galaxies may be described by an exponential timescale of ~ 6 Gyrs up to z ~ 1.5, and proceeds much faster 
thereafter, with an exponential timescale of ~ 0.6 Gyrs. We also show that there is some evidence for a differential evolution of 
the Galaxy Stellar Mass Function, with low mass galaxies evolving faster than more massive ones up to z ~ 1 — 1.5 and that 
the Galaxy Stellar Mass Function remains remarkably flat (i.e. with a slope close to the local one) up to z ~ 1 — 1.3. 
Conclusions. The observed behaviour of the Galaxy Stellar Mass Function is consistent with a scenario where about 50% of 
present-day massive galaxies formed at a vigorous rate in the epoch between redshift 4 and 1.5, followed by a milder evolution 
5_l ' until the present-day epoch. 
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1. Introduction 

The observational evidence of the continuous increase 
of the stellar content of galaxies over cosmic times has 
emerged only recently from a series of observations 
and surveys which have made use of new sensitive IR 
instrumentation. Following the first pioneering stud- 
ies HGiallongo et al. T998| IBrinchmann fc Ellis 20001 
|Papovicn^t™al~20CTl| that set-up the technique for 
estimating stellar masses in high redshift galaxies, several 
surveys pointed out that the global stellar content of the 
Universe, as measured by the average stellar mass density, 
increases with cosmic time ( Brinch mann &: Ellis 2 000 
Dickinson et al. 20031 IFontana et"al~2 003 
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Rudnick et al. 20031 IGlazebrook et al. 20041 

Drory et al. 2004| IFontana et al. 20041 F04 hereafter, 
Rudnick et al. 2006|> . 

While the average stellar mass density provides a 
global picture of the process of stellar assembly, the 
Galaxy Stellar Mass Function (GSMF in the following) 
provides a more detailed view on how this process evolves 
as a function of the galaxy mass itself. At low red- 
shift, accurate GSMF have been obtained from the 2dF 
l|Cole et al. 2001 |l or 2MASS-SDSS <|Bell et a I. 2003II sur- 
veys. Because of the large area and of the accuracy in 
the redshift estimate that is needed to compile a GSMF, 
very fe w GSMFs have bee n obtained so far. The MUNICS 
survey ( |Drory et al. 2004| ) and the K20 survey (F04) first 
explored the evolution of the GSMF up to z ~ 1 — 2 in 
a somewhat complementary fashion. The former adopted 
a wide, relatively deep sample on a wide area (~ 5000 
objects at K < 19), mostly relying on photometric red- 
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shifts, sampling the GSMF up to z ~ 1.5, while the latter 
adopted a deeper but smaller sample (~ 500 objects at 
K < 20), with excellent spectroscopic coverage, extending 
up to z = 2. In the overlapping redshift ranges, the two 
GSMF agree quite well, and both analyses suggested a de- 
cline in the density of massive galaxies at z ~ 1, of the 
order of 50-70% with respect to the local value. Recently, 
|Drory et al. 2005| extended such analysis to higher [z ~ 5) 
redshift, adopting a combination of / and -fT-selected sam- 
ples. 

F04 also pointed out a tentative evidence of a differ- 
ential evolution of the GSMF, such that massive galaxies 
appear to evolve less than low mass galaxies, at least up 
to z ~ 1.5. As we shall show in the following, this trend 
is also shown by our new data. 

Along a different line, evidence has emerged that the 
GSMF for galaxies of different spectral or morphological 
types shows a significant evolution, with an increase in 
the fraction of stellar mass residing in late type galax- 
ies at z ~ 1 — 1.5 with respect to low redshift (F04, 
|Bundy et al. 2005| IFranceschini et al. 2006|> . 

Although these surveys have already provided a first 
picture of the evolution of the GSMF, the uncertainties 
involved in this exercise are still large. In particular, the 
exploration of the high redshift Universe is still very lim- 
ited. In several cases, the lack of long-wavelength data 
hampered the estimate of the stellar mass at high red- 
shift, for which an adequate sampling of the rest frame 
optical-near infrared part of the spectrum is essential. In 
addition, the collection of a large sample of high redshift 
galaxies requires a combination of large areas and deep 
near-IR observations. For different reasons, all the sur- 
veys mentioned above suffer of these limitations, although 
to a different extent. 

In this context, the GOODS-South survey provides an 
excellent opportunity of improving in a significant way 
our knowledge of the high-z GSMF. The combination 
of deep, wide IR observations with Spitzer (in its four 
channels from 3.5 to 8 /zm) and VLT-ISAAC (in J, H 
and Ks), coupled with high quality imaging data in the 
optical domain (both ACS and VLT-VIMOS), and of a 
large, extended spectroscopic coverage make the GOODS- 
South field ideal for this investigation. From this pub- 
lic data set we have obtained a multicolour catalog of 
faint galaxies, that we named GOODS-MUSIC (GOODS 
MUltiwavelength Southern Infrared Catalog), that we de- 
scribe in G razian et al. 2006al 

In this work we use the GOODS-MUSIC sample to 
improve the previous estimates of the GSMF in several 
ways. The most important is that our analysis includes 
the 3.5 — 8i/m Spitzer observations of the complete Ks- 
selected data set. In addition, it extends to Ks magnitudes 
deeper than any previous survey, enabling us to obtain a 
complete sample of galaxies up to z ~ 4, and to measure 
the slope of the low mass side of the Galaxy Stellar Mass 
Function up to z ~ 1.3. Compared to our previous analysis 
of the K20 survey data, the present data set provides a 
final sample that is 6x larger than the K20 one, on which 



we adopt a more sophisticated technique to parametrize 
the evolution of the Galaxy Stellar Mass Function with a 
redshift-dependent Schechter fit, that provides interesting 
clues on the evolution of massive galaxies. 

The paper is organized as follows. In Sect. 2, after re- 
minding the basic feature of our dataset and of the proce- 
dure that we adopt to extract the stellar masses from each 
galaxy, we discuss how the inclusion of the Spitzer bands 
affects the mass estimates and highlight the selection ef- 
fects that will affect our analysis. In Sect. 3, we present the 
basic results of our analysis, namely the Galaxy Stellar 
Mass Function and the resulting mass density, both in 
a binned as well as in a parametric fashion. In Sect. 4, 
we compare our basic findings with the prediction of re- 
cent theoretical models. Finally, in Sect. 5, we focus on the 
highest redshift range, to discuss the reliability of the pho- 
tometric redshifts on a few, intriguing objects that might 
be at very high redshift. The discussion and conclusions 
are summarized in Sect. 6. 

All the magnitudes cited below are in the AB sys- 
tem. To scale luminosities and compute volumes we have 
adopted a "concordance" cosmological scenario with = 
0.7 and H = 70 km s^ 1 Mpc" 1 . 

2. Stellar Masses in the GOODS-MUSIC sample 

2.1. The Data 

We use GOODS-MUSIC, a multicolour catalog extracted 
from the deep and wide survey conducted over the 
Chandra Deep Field South (CDFS) in the framework 
of the GOODS public survey. The procedures that we 
adopted are described at length in IGrazian et al. 2006a 
We remind here only the basic features. 

The data comprise a combination of images that ex- 
tend over 14 bands, namely U-band data from the 2.2ESO 
(U 35 and U 38 ) and VLT-VIMOS {U V imos), the F435W, 
.F606W, F775W and F850LP (Z 850 ) ACS images, the 
JHKs VLT data and the Spitzer data provided by IRAC 
instrument (3.6, 4.5, 5.8 and 8.0 jim). From this dataset, 
we have obtained a multiwavelength catalog of 14847 ob- 
jects, selected either in the ifgso and/or in the Ks band. 
For the purposes of the present work, we will mainly use 
the i^s-selected catalog. This consists of 2931 galaxies 
(after removal of known or candidate AGNs and Galactic 
stars), 1922 of which have Uvimos coverage, 1762 have H 
coverage, and all have a complete coverage in the remain- 
ing 12 bands, most notably including the Spitzer ones. 
Since the detection mosaics have an inhomogeneous depth, 
we have divided the Ks sample into 6 independent cata- 
logs, each with a well defined magnitude limit and area, 
that we use to compute mass functions and mass densities. 
However, the largest fraction of the sample has a typical 
magnitude limit of Ks ~ 23.5 that we will adopt for more 
qualitative arguments. 

Colours have been measured using a specific soft- 
ware for the accurate "PSF-matching" of space and 
ground based images of different resolution and depth, 
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that we have named ConvPhot (De Santis et al 2006). 
We have cross correlated our catalog with the whole 
spectroscopic catalogs available to date, from a list 
of surveys, assigning a spectroscopic redshift to more 
than 1000 sources. We note that in this work we use a 
spectroscopic sample that is wider than that presented 
in IGrazian et al. 2006al thanks to the increased num- 
ber of spectra publicly available QVanzella et al. 2*0 06). 
Finally, we have applied our photometric redshift 
code, developed and tested over the years in a series 
of works IjFontana et al. 20001 ICimatti et al. 20021 
IFontana et al. 20031 IFontana et al. 20041 

|Giallongo et al. 2 005), that adopts a standard x 2 
minimization over a large set of templates obtained 
from synthetic spectral models. The comparison 
with the spectroscopic sample QGrazian et al. 20 06a 
IGrazian et al. 20"0 6b) shows that the quality of the 
resulting photometric redshifts is excellent, with a r.m.s. 
scatter in Az/(1 + z) of 0.03 and no systematic offset. 

In summary, the final sample adopted in this work con- 
sists of 2931 galaxies, complete to a typical magnitude of 
Ks ~ 23.5, over an area of 143.2 sq. arcmin, 815 of which 
with reliable spectroscopic redshift (i.e. 28% of the total 
sample) and the remaining fraction with well tested 14 
bands photometric redshifts. 

2.2. Stellar Masses in the Spitzer era 

The method that we applied to estimate the galaxy stellar 
masses on this dataset is exactly the same that we devel- 
oped in previous papers fFo ntana et al. 20031 F04), and 
similar to those adopted by other groups in the literature 
(e.g. IDickinson et al. 20031 |Drory et al. 2004| ). Briefly, it 
is based on a set of templates, computed with standard 
spectral synthesis models (Bruzual & Chariot 2003 in 
our case), chosen to broadly encompass the variety of 
star-formation histories, metallicities and extinction of 
real galaxies. To compare with previous works, we have 
used the Salpeter IMF, ranging over a set of metallic- 
ities (from Z = O.O2Z to Z = 2.5Z Q ) and dust ex- 
tinction (0 < E(B — V) < 1.1, with a Calzetti extinc- 
tion curve). Details are given in Table 1 of F04. For each 
model of this grid, we have computed the expected magni- 
tudes in our filter set, and found the best-fitting template 
with a standard x 2 normalization. The stellar mass and 
other best-fit parameters of the galaxy are found after 
scaling to the actual luminosity of the observed galaxy. 
Pros and cons of the method have been discussed in 
several papers UPapovich et al. 2001| IFontana et al. 20041 
|Shapley et al. 2004| ), and we refer to them for a detailed 
discussion of the systematics involved in this exercise. 

The major difference with all previous estimates of the 
GSMF arises from the inclusion of the 4 Spitzer bands, 
longward of 2.2 /j,m. For galaxies at z > 2, these bands are 
essential to sample the spectral distribution in the rest- 
frame optical and near-IR bands, that are necessary to 
provide reliable constraints on the stellar mass. We quan- 
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Fig. 1. Relation between the stellar masses derived with 
the inclusion of the Spitzer bands (14 bands) and those 
without the Spitzer bands (10 bands), as a function of red- 
shift, on the Ks-selected galaxies of the GOODS-MUSIC 
sample. In the lower panel we plot the ratio between the 
two estimates, for each galaxy in the -ftTs-selected sample. 
In the upper, the distribution of such a ratio in different 
redshift ranges. 



titatively assess their importance in Fig^ where we plot 
the ratio between the mass estimates with (14 bands) and 
without (10 bands) the Spitzer bands on the sample of 
ifs-selected galaxies, as function of redshift. 

It is immediately appreciated that the inclusion of 
the Spitzer bands does not provide statistically significant 
changes at z < 2, where most of the galaxies are fitted with 
the same best-fitting models in both cases. The scatter in 
the mass estimate is entirely consistent with the expected 
uncertainty due to model degeneracy. At higher redshift, 
the scatter increases significantly. At z ~ 2.5, the r.m.s. 
fluctuation is of about 0.23, although with no systematic 
shift. At higher redshift, the r.m.s. is as large as 0.45 dex, 
with evidence of systematic overestimate (of 0.2 dex) of 
the stellar mass when the Spitzer bands are not included. 

The same improvement can be found by looking at the 
uncertainty in the mass estimate. As in F04, we computed 
the la confidence level on each mass estimate by scanning 
the x 2 levels, allowing the redshift to change in case of 
objects with photometric redshift. In Fig(2] we compare 
these estimates when the Spitzer bands are included (14 
bands) or not (10 bands). As clearly shown, the formal 
intrinsic uncertainty decreases from a typical value of 60% 
to 40% on the global sample, and decreases by a factor of 
three for objects at z > 2. 
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Fig. 2. Effect of the Spitzer bands on the mass estimate. 
Upper panel we plot the distribution of the relative un- 
certainty (at the la level) in the mass estimate, when the 
Spitzer bands are included (14 bands) or not (10 bands). 
The distribution is computed over all i-Ts-selected galax- 
ies. Lower panel: Ratio of the uncertainty (at the la level) 
in the mass estimate, when the Spitzer bands are included 
(14 bands) or not (10 bands), as a function of redshift. 



These results confirm that, in the absence of the 
Spitzer bands, the estimates of the stellar mass for galaxies 
at z > 3 are very uncertain and possibly biased, such that 
detailed astrophysical analysis based on such estimates are 
likely premature. 

2.3. Stellar masses of the GOODS-MUSIC data sample 

In Fig|3| we plot the resulting stellar masses of the 
GOODS-MUSIC sample, as a function of redshift. We 
include in the plot both the objects of the if s-selected 
sample as well as those of the Zgsn-selected one that have 
detected flux in the K s band. 

One can immediately appreciate that massive galax- 
ies (i.e. those above the z = characteristic mass M* ~ 
10 11 M Q ) are thoroughly detected up to at least z ~ 4, 
and possibly even at higher z, although the reliability of 
photometric redshifts for the massive objects at z ~ 5 is 
uncertain, as discussed in Section 5. We also note that such 
plot appears to be qualitatively different from the analo- 
gous plot by |Drory et al. 20 05 that contains several quite 
massive galaxies (close to 10 r2 M & ) up to z ~ 5. This dif- 
ference is likely due to a combination of effects, induced by 
the different techniques applied for photometry and pho- 
tometric redshifts, the availability of an extended spectro- 
scopic coverage, and probably even more important the 
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Fig. 3. Stellar Masses for individual galaxies of the 
GOODS-MUSIC sample, as a function of redshift. Upper 
panel: the ifs-selected (filled circles) and the Z S50 - 
selected sample (small dots) on the overall redshift range. 
Lower panel: same quantities, for galaxies at z > 2. 
Objects in the -ftfs-selected sample are shown with large 
empty circles if they are also included in the Z 8 5o-selected 
one, and with large filled circles if included only in the 
ifs-band catalog. Small dots are from the remaining Z^q— 
selected objects. In both panels, lines represent the com- 
pleteness threshold, computed as described in the text. 
Thick lines correspond to the completeness threshold for 
a ifs-selected sample with respect to passively evolv- 
ing models (continuous line) and dusty starburst with 
E(B — V) = 1.1 (dashed). The thin lines show the same 
thresholds if a 4.5^im-selected sample would be used, for 
the two same spectral models. 



lack of Spitzer observations in the |Drory et al. 2005| sam- 
ple. 

In the lower panel of FigOU we focus on the high red- 
shift range of our sample, and we use different symbols 
to differentiate among galaxies that are selected only in 
the K s-selected sample, those that are selected both in 
Ks and .Zsso an d those that are selected in Z &50 only. 
As expected, galaxies that are selected in Z^q only, and 
hence have a very low flux in Ks, populate the low mass 
region of the distribution. On the contrary, a number of 
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galaxies selected only in Ks contribute to the population 
of massive galaxies at z > 2: these are typically very faint 
or sometimes even undetected in the Z S50 -band, resulting 
from highly dust-reddened or passively fading stellar pop- 
ulation at high redshift. These objects significantly con- 
tribute to the number density of massive systems. In a 
similar way, one expects that other massive galaxies at 
these redshifts might be included by selecting the sample 
at even higher wavelengths, as those obtained with the 
Spitzer observations used here. 

2.4. Incompleteness and selection criteria 

It is important to remark that our sample, as well as any 
other magnitude -limited sample, does not have a well 
defined, sharp limit in stellar mass. This incompleteness 
effect arise from the fact that galaxies have a range of M/L 
ratios, and can be illustrated as follows. In the GOODS- 
MUSIC sample, at z — 1 we find that galaxies have a 
range of M/L extending from 0.9 (for redder objects) to 
0.046 (for bluer objects). Assuming for simplicity a sharp 
magnitude limit Ks — 23.5, this corresponds to stellar 
masses from M@ ~ 1.5 • 10 10 (for redder objects) to M Q ~ 
6 • 10 8 (for bluer objects). However, since galaxies have a 
range of M/ L ratios, additional objects with stellar masses 
in the range 6 • 10 8 < M Q < 1.5 • 10 10 will lay at fainter 
fluxes, Ks > 23.5, which are not included in the present 
sample. Hence, the census of galaxies in the range 6 • 10 8 < 
Mq < 1.5-10 10 will definitely be incomplete in our sample, 
at z ~ 1. A complete visualization of this effect can be 
obtained by looking at Fig 16 of F04. Needless to say, this 
effect exists in any magnitude-limited sample and at any 
redshift, although going to redder wavelengths alleviates 
its impact. 

This effect has been described at length by Dickinson 
et al 2003, Fontana et al 2003 and F04, and can be dealt 
with in two different ways. First, one can compute a lim- 
iting threshold in mass, such that one can obtain a well 
defined mass-selected sample. All the relevant statistics 
(mass densities, mass functions etc) should be obtained 
using only the objects above such limit in stellar mass. 

To better exploit the statistics of the sample, we have 
developed a technique to correct for the incomplete cover- 
age of M/L ratios, that allows to at least partially recover 
the fainter side of the sample. This technique will be used 
to estimate the GSMF and will be described in Sect. 3. 
In our sample, this technique allows to extend the GSMF 
by about OAdex in each redshift bin: the leading effect 
remains therefore the conservative estimate of the com- 
pleteness limit, that we have derived as follows. 

The limiting threshold in mass can be obtained by 
computing the maximal M/L that is allowable at each 
redshift and multiplying it by the magnitude limit of the 
survey. Below such a threshold in mass, one can still de- 
tect galaxies and measure their mass, but will start to miss 
galaxies of the same mass but with lower luminosity (i.e. 
higher M/L). Unfortunately, such a limiting threshold de- 



pends on the assumed spectral properties for the targeted 
galaxies (and obviously on the bandpass adopted to select 
the catalog). 

In the simplest case, one can compute the maximal 
M/L for passively evolving systems. In this case, it is easy 
to draw a selection curve, that we plot in Figdl (thick solid 
line), that we computed with a maximally old, single burst 
model normalized to Ks = 23.5. However, large values of 
M/L can also be found in heavily extincted star-forming 
galaxies. To estimate this effect, we have adopted a simple 
star-forming model with a variable amount of extinction, 
adopting a Calzetti extinction law. It turns out that, for 
E(B — V) ~ 0.5, the corresponding selection curve roughly 
corresponds to the "passively evolving" one, and we do not 
plot it for simplicity. Increasing the amount of dust, the 
selection curve shifts to higher masses. In Figdlwe plot the 
case for E{B — V) ~ 1.1 (thick dashed line), that is among 
the highest observed in spectroscopically confirmed EROs 
l|Cimatti et al. 2002jl . 

The result must be read as follows: our ifs-selected 
sample is expected to be complete against passively evolv- 
ing objects with mass above the thick solid line in Fig|3 
i.e. grossly with mass M* > 10 11 M Q up to z ~ 4; however, 
at this mass level it becomes progressively biased against 
the detection of extremely dusty objects, probably already 
at z > 2. 

It is interesting to predict what would be the advan- 
tages of adopting a Spitzer-selected sample, scaled to the 
expected sensitivity of the GOODS survey. We have com- 
puted the same selection curves (for passively evolving 
and dusty objects) adopting a magnitude limit at 4.5/im 
of m4.5 < 23.5, and we have overplotted them in FigO 
as thin solid and dashed lines. A first result is that the 
adoption of a 4.5/Lim-selected sample would not allow to 
extend the sampling of passively evolving galaxies at lower 
masses than in our catalog, but rather to extend the se- 
lection of objects at the same level of M* ~ 10 11 M Q well 
beyond z ~ 4. Probably the most important effect is that 
the adoption of a 4.5/mi-selected sample would allow the 
detection of heavily extincted, dusty objects well beyond 
z ~ 2, probably up to z ~ 3 at the M* ~ lO n M0 level. 
Objects detected in the 4.5/im band, but very faint or even 
undetected in Ks, do exist in the GOODS area (see Yan 
et al 2005 for a few cases): unfortunately, their inclusion 
in the present analysis would require a detailed analysis 
(most notably of the reliability of their photometric red- 
shift), that would go beyond the purposes of this work and 
that we defer to a dedicated paper. 



3. The Galaxy Stellar Mass Function 

3.1. Computing the Mass Function 

We have used the ifs-selected sample described above to 
compute the GSMF at various redshifts. The GSMF is 
computed both by using the standard l/V max formalism 
as well as by fitting a Schechter function to the unbinned 
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Fig. 4. Galaxy stellar Mass Functions in the GOODS-MUSIC sample, in different redshift ranges. Big circles repre- 
sent the Galaxy Stellar Mass Functions of the ifs-selected sample, computed with the 1/V ma x formalism up to the 
appropriate completeness level, as described in the text, while small triangles show the Galaxy Stellar Mass Functions 
of the ^so-selected sample, computed without any correction for incompleteness coverage of the M/L at faint limits. 
The dashed region represents the local GSMF of lCole et al. 20011 encompassing its original and the "rescaled" version 
that we obtained in F04 (see text for details). The solid line is the evolutionary STY fit of our data, computed over 
the global redshift range 0.4 < z < 4, that we describe in Section 3.2. The observed mass functions are available in 
electronic form at http://lbc.oa-roma.inaf.it/goods/massfunction 



data. The results of the former method are described in 
Sect. 3.2, while the latter in 3.3. 

The GSMF must be computed taking into account the 
incompleteness effects described above. To correct for this 
effect, we have introduced a correction technique, which 



is carefully described in F04, and that we remind briefly 
here. We start from the threshold computed from pas- 
sively evolving system, i.e. the thick solid line of FigOU 
below such completeness threshold only a fraction (/obs) 
of objects of given mass will be observed. We then obtain 
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(at any redshift) the observed distribution of M/L ratios 
for objects close to the magnitude limit of the sample. 
Then, using such distribution, we compute for each mass 
and redshift the fraction 1 — / b s of galaxies lost by ef- 
fect of the incomplete coverage of the M„/Lk ratio, by 
which we obtain Jobs- As shown in F04, Appendix B, rela- 
tively simple analytic formulae can be used to describe this 
process. We adopted the same technique to the GOODS- 
MUSIC sample, verifying that these simple analytic ex- 
pression still provide a good fit to the observed M*/Lk 
distribution, despite the greater depth and area of the 
GOODS-MUSIC sample. The correction is applied only 
for a limited range of masses, where the correction is less 
than 50% (i.e. / b s > 0.67. This correction factor is then 
applied to the volume element V max of any galaxy in the 
1/Vniax binned GSMF as well as in the number of detected 
galaxies that enters in the Maximum Likelihood analysis 
used to obtain the Schechter best-fits. As shown in Fig. 17 
of F04, the inclusion of this treatment allows to extend the 
computation of the GSMF to somewhat (OAdex) smaller 
masses than by adopting the strict completeness limit. In 
contrast, neglecting the whole treatment of incomplete- 
ness introduces a remarkable and unphysical bending of 
the GSMF at low masses. 

A further somewhat technical point is related to the 
choice of the local GSMF that is used to estimate evolu- 
tion. As described in F04, the parameter grid used as input 
to the spectral synthesis model adopted here is different 
from the one adopted in the derivation of the local GSMF 
by ICole et al. 20011 in particular for the absence of any 
constraint on the galaxy ages. This systematic difference 
may lead to a small but systematic bias in the measure of 
galaxy masses, that can be described by an average rela- 
tion (log(M) = 1.027 * log( M co ; e ) - 0.3955). By applying 
this relation to the original ICole et al. 20011 local GSMF 
we obtain a "rescaled" version of the local GSMF, where 
differences are in practice noticeable only on the exponen- 
tial tail of the GSMF. For sake of clarity, we will in the 
following represent both the original Co le et al. 2001l local 
GSMF as well as its rescaled form that we derived in F04. 

3.2. The distribution of massive galaxies up to z~ 4 

The Galaxy Stellar Mass Function obtained in the 
GOODS-MUSIC sample is shown in FiglU where it is 
computed in several redshift bins. We also plot in FigQ] 
the GSMF obtained on the .Zsso -selected sample, naively 
computed with no correction for the incomplete cover- 
age of the M/L distribution. We explicitly note that such 
a Zg5Q-selected GSMF is quite consistent with the Ks- 
selected one at low and intermediate redshifts, before fad- 
ing at the lowest mass end because of the lack of correc- 
tion for incompleteness. Conversely, it misses a significant 
fraction of massive galaxies at z > 2 — 3, as expected in a 
Zgso-selected sample. 

We also remark that the error bars of Fig0]have been 
computed with a full Monte Carlo simulation where we 



take into account the redshift probability distribution of 
each galaxy in the sample. 

Three different results can be inferred even from a 
visual inspection of FigQ] First, the density of massive 
galaxies clearly decreases with redshift, relatively mildly 
up to z ~ 1.5 and then much more convincingly at z > 2. 
At z ~ 3.5, their density is at least a factor 10 x lower 
than in the local Universe. This result confirm the trend 
that we already noted in F04 (see also Drory et al 2005). 

Second, the number density of lower mass system 
evolves significantly with redshift (at z ^ 1 the density 
of galaxies with M ~ 1O 1O M0 is 4 times lower than the 
local one), while the slope of the GSMF is relatively flat 
up to at least z ~ 1 — 1.5, and does not appear to steepen 
significantly with respect to the local one. The latter result 
is apparently at variance with the results of in F04, where 
we found a significant steepening of the GSMF at z ~ 0.7, 
with the slope index a changing from a — —1.18 at z = 
to a = —1.4 at z = 0.7. We first carefully checked that no 
systematic effect is responsible for this result: in particu- 
lar, we verified that the GSMF computed with our new 
sample in the very area of the K20 survey is comparable 
to our early results (actually it is nearly perfectly coinci- 
dent), despite the new photometry of the whole sample. 
We furthermore note that even in the GOODS-MUSIC 
sample the slope of the GSMF is steeper than the average 
in the redshift bin z = 0.8 — 1, as also shown by the para- 
metric analysis presented in Sect. 3. 3, probably due to a 
cosmic variance effect. Taking also into account the agree- 
ment with the Z 85 o-selected GSMF, we conclude that our 
new sample, that is much deeper and wider than the K20 
sample, allows to better represent the global evolution of 
the slope of the GSMF with redshift. 

Finally, there is some evidence of a differential evolu- 
tion of the GSMF, in the sense that more massive galax- 
ies appear to evolve less than low mass galaxies, at least 
up to z ~ 1.5. This trend was already noted in F04, but 
the more limited statistics prevented a robust quantitative 
conclusion. 

Although these conclusions are clearly supported by 
our data, we are aware that the 1/V max approach is very 
sensitive to the biases induced by Large Scale Structures 
(LSS) and clustering, and that our sample is clearly af- 
fected by these effects. A well known effect is the ap- 
parent steepening of the faintest points in case of clus- 
tering ( |Heyl et al. 1997j ): this is present in our sample at 
1 < z < 1.3 and 1.6 < z < 2, due to two large over- 
densities at the lower limit of the bin. Such effects are 
amplified by the arbitrariness in the choice of the red- 
shift bins, and by our choice of keeping a small binning 
factor in mass. If one compares our sample with the lo- 
cal GSMF of ICole et al. 20011 it is possible to see that 
fluctuations in the number of massive galaxies exist, su- 
perimposed to the global evolution of the GSMF. In par- 
ticular, the number of massive galaxies and even the de- 
tection of galaxies with M* > 10 n M Q is higher in the 
redshift bins where LSS have been detected in the CDFS 
field, e.g. at z = 0.67,0.73, z ~ 1.61 l|Gilli et al. 20031 
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IVanzella et al. 20051 IVanzella et al. 2006jl . while our sam- 
ple is relatively devoid of massive galaxies at z ~ 0.9, 
where an underdensity of galaxies exists in the CDFS 
IjPerez-Gonzales et al 200 5). To overcome these limita- 
tions, we have developed the STY analysis that we shall 
describe in the following section. 

3.3. The STY approach 

An independent approach to the evaluation of the GSMF 
is the STY fitting method ( |Sandage et al 1979] ), mutil- 
ated from the corresponding formalism of the Luminosity 
Function. This method assumes an analytic expression 
of the GSMF and derives the free parameters with a 
Maximum Likelihood analysis. It is less sensitive to clus- 
tering and LSS and provides more quantitative hints on 
the evolution of the global population. To minimize the 
impact of LSS, and to avoid any arbitrariness in the defi- 
nition of the redshift bins, we use here a global, redshift- 
dependent parametrization of the GSMF. At each redshift, 
the number density <fi(M, z) of galaxies with mass M* is 
assumed to have the functional form of a Schechter func- 
tion, 

4? (z) log (10) [10< M - M * W) exp(-10( M - M * W )) (1) 

where M — logio(M*/Mo) and the free parameters 
4>*{z), a*(z) and M*(z) are functions of redshift. For these 
parameters, we have found that the following simple rela- 
tions provide an adequate fit to the overall evolution, in 
the redshift range 0.4 < z < 4: 



a (z) = a Q + a 1 ■ z 

M*(z) = M* + M* ■ z + M* ■ z 1 



(2) 
(3) 
(4) 
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Fig. 5. Contour levels of all the pairs of free parameters 
involved in the Schechter fitting of the GSMF, as obtained 
by a Maximum Likelihood minimization. Inner ellipses 
represent the la contour, outer ones the 3<r contour. All 
contours have been obtained with marginalization of the 
other free parameters. 



A first robust result of the STY analysis is that the 
slope of GSMF is remarkably flat and with a small redshift 
evolution: up to z ~ 1 — 1.3, the highest redshift where 
the slope is reliably measured in our sample, the slope a 
changes of about 0.1 in our fiducial best-fit model. As dis- 
cussed before, this is accompanied by a sensible decrease 
in the number density of low mass galaxies, that for galax- 
ies with M ~ 10 w M(t, is 4 times lower than the local one. 



where </>q, <p\ , aj, al , Mq, M*, M| are free parameters. 

In this formalism, the three zero-th order parameters 
(Mq , (j)Q and «g) should ideally reproduce the local values, 
as estimated for instance bv lCole et al. 200ll For this rea- 
son we fixed the local characteristic mass and slope to the 
local values of lCole et al. 200Tl and derived the other pa- 
rameters. The resulting fit is displayed in FigQ]as a solid 
line, and the relevant parameters are listed in Table ^ 
The Maximum Likelihood analysis allows also to perform 
an estimate of the error budget on the fitted parameters: 
we show in FigEl the la and 3cr contour levels on each 
pair of free parameters, obtained after marginalization of 
the other parameters. It is clearly shown that all param- 
eters are reasonably well constrained, in particular those 
(M*,M% and a*) that are more important for the phys- 
ical conclusions that we will draw in the following. The 
resulting evolution with redshift of the characteristic mass 
M*(z) and of the slope a*(z) is shown in FigEl where they 
are also compared with the corresponding values obtained 
by a Schechter fit within each redshift bin. 



The evolution of the characteristic mass M*(z) is more 
complex. According to the best-fit values, the evolving 
M*(z) is fitted with a second order law with negative cur- 
vature (M 2 * = -0.07±0.02) and a peak at about z ~ 1.35, 
such that the characteristic mass of the GSMF initially 
increases, reaches a maximum around z ~ 1.35 and than 
decreases. 

We carefully verified that this behaviour is not a spu- 
rious result of our choice of the functional form of Eq|3] 
where the second order law for M*(z) might introduce 
such a non-monotonic behaviour. First, as shown in Fig0 
the second order term M 2 is definitely negative, suggest- 
ing that M*(z) actually increases up to z ~ 1 — 1.3 and 
decreases thereafter. Such a behaviour is also consistent 
with the trend observed in the values of M* in the indi- 
vidual redshift bins, as shown in FigU3 More important, 
we arrived at this choice after discarding simpler solutions. 
In particular, we tested a linear form for M*(z) (i.e. with 
Ml = 0), as well as a more traditional logarithmic func- 
tion (M*(z) = Mq + M* ■ log(l + z)). The output of the 
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Redshift 

Fig. 6. The evolution of the Schechter parameters M* 
(characteristic mass) and a (slope) as a function of red- 
shift. Points represent the individual fits to the single red- 
shift bin of Fig0] with their uncertainty, while lines rep- 
resents the evolutionary STY fit of our data, computed 
over the global redshift range 0.4 < z < 4, that we de- 
scribe in Section 3.2. Solid line is the second-order form 
for M* (our fiducial model); dashed line is for M* = 
Ml x z + M| x z M 3 ; dotted line is for M* oc log(l + z). 
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Fig. 7. Ratio between the observed Galaxy Stellar Mass 
Functions in the GOODS-MUSIC sample and the local 
GSMF, in the same redshift bins of FigQ] The observed 
data are normalized to two different choices of the local 
GSMF, namely the original one of ICole et al. 20011 (filled 
squares) and its rescaled version, that we described in the 
text (filled circles). The shaded area shows the correspond- 
ing ratios for the global Schechter fit, that we describe in 
Sect. 3.2. 



latter form is also shown in Fig|5] In all these cases, we 
were not able to fit the massive side of the GSMF at inter- 
mediate redshifts (0.6 < z < 1.5), that is systematically 
undcrfitted, by nearly 50%. As a further test, we have also 
allowed for the exponent of the high order term to be free 
(i.e. M*(z) = Mq + Ml ■ z + M| • z M 5), and found that it 
is definitely larger than 1 (M| = 1.2 ±0.05), providing an 
evolution that is similar to our fiducial model, as shown 
by the dashed line in FigHj] We conclude that such evolu- 
tionary form is a robust property of our sample, and we 
shall adopt the second order form for M*(z) shown in Eq. 
4 as our best-fit fiducial model. 

Other warnings arise from possible biases in our sam- 
ple. The most obvious sources of biases are sample vari- 
ance and the use of photometric redshifts for a relatively 
large fraction of the sample. We expect the latter uncer- 
tainty to be minor, given the large spectroscopic coverage 
of the GOODS survey for the brightest galaxies (where 
most of the evolution of M* is measured) and the good ac- 
curacy that we achieve in photometric redshifts. Sample 
variance is a more serious concern, since the GOODS- 
South field is definitely too small to avoid spurious effects 
due to over- and under-densities along the line of sight. 
In the redshift interval 0.6 < z < 1 a range of differ- 
ent densities exist, from overdensities (around z = 0.7) 
to underdensities (at z ~ 0.9), and the increase of M* is 



reassuringly seen in both. However, much wider surveys 
with comparable depth are definitely required to confirm 
or dispute our finding. 

Turning to physical interpretation of this result, we 
remark that the evolution of the characteristic mass M* 
reflects an evolution of the shape of the GSMF, that is 
progressively skewed toward larger masses as it evolves 
from z = to z ~ 1.3. The parametric analysis of this 
section enforces what we already highlighted in Sect. 3. 2, 
i.e. that the GSMF appear to have a differential evolution, 
with more massive galaxies evolving less than low mass 
galaxies, up to z ~ 1 — 1.3. We note that this does not 
imply that the mass density increases, since the increase 
of M*(z) is counterbalanced by a decrease of the overall 
normalization <fi*(z). 

The issue of the differential evolution of the GSMF 
is potentially very important, since it is directly related 
to the "downsizing" scenario for galaxy evolution, such 
that it is worth a more careful examination. We plot in 
Fig'0the ratio between our GSMF and the local one, for 
different redshifts. In addition to the obvious overall evo- 
lution of the GSMF with redshift, it is shown that the 
density of galaxies above 1O 11 M - albeit decreasing - 
remains closer to the local value up to z ~ 1.15 than 
that of the lower mass population. Such a trend appears 
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Table 1. Best fit Schechter parameters 
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to be progressively stronger for the more massive galax- 
ies of our survey, i.e. those with a stellar mass in excess 
of 3 x 1O 11 M . At larger redshifts, however, this trend 
eventually breaks down, such that the density of massive 
galaxies undergoes a strong evolution in the redshift range 
1.5 < z < 3. 

3.4. The integrated Stellar Mass Density 

A more global view on the rise of the galaxy mass as a 
function of cosmic time is provided by the Stellar Mass 
Density (SMD), that is obtained by integrating the GSMF 
over all masses (we choose in particular to integrate from 
log(M*/M ) = 8 to log(M*/M ) = 13). This is displayed 
in Fig|Sl where we plot as a solid line the SMD as ob- 
tained by integrating our Schechter fiducial model. At 
z > 1.3, where we do not adequately sample the faint 
end of GSMF, we rely on the extrapolation of the slope 
a provided by our fiducial model, that corresponds to a 
very mild steepening with redshift. In Tab[21we provide 
the SMD values as obtained by integrating the Schechter 
fits as well as those directly observed. When compared 
also with other surveys, FiglHldepicts a scenario where the 
global stellar mass density has evolved relatively slowly 
(i.e. by a factor of two) over the last 8 Gyrs, and more 
rapidly at higher redshift, with a decrease by an order of 
magnitude at z ~ 3. 

A more direct result of our survey is the SMD in mas- 
sive galaxies, defined as those above M = 10 11 M Q , that 
we directly detect in our survey up to z ~ 4. 

As we show in the lower panel of Fig|SJ the evolution 
of the stellar mass density of massive galaxies increases 
fastly over the first 3-5 Gyrs in the history of Universe, 
and thereafter proceeds at a much slower pace. It is possi- 
ble to grossly reproduce this behaviour by an exponential 
law, p(z) oc e~*/ r , characterized by two different e- folding 
timescales. At high redshift (z > 1.5, i.e. look back time 
> 9 Gyrs), the timescale that we derive from our data is of 
the order of 0.6 Gyrs. At later cosmic times, the timescale 
is at least a factor of 10 larger: we obtain 6 Gyrs from 
z = to z = 1.5. 




5 10 
Look Back Time (Gyr) 

Fig. 8. Global Stellar Mass Density as function of red- 
shift. Upper panel: the total Stellar Mass Density (inte- 
grated from log(M*/M©) = 8 to log (M»/M Q ) = 13) as es- 
timated in the present work, obtained by an integration of 
our fiducial Schechter function (solid line). It is compared 
with other estimates adopted from the literature: the local 
value by Cole et al. 2001 (large empty circle), Fontana et 
al. 2004 (red filled circles), Fontana et al. 2003 (HDFS, 
open diamonds), Dickinson et al. 2003 (open squares), 
Cohen et al. 2001 (open crosses), Brinchmann and Ellis 
2000 (filled diamonds), Drory et al. 2005 (open triangles, 
K-selected sample with horizontal basis, I-selected sam- 
ple with vertical basis). To avoid further confusion, we 
do not plot here the mass density computed in this pa- 
per. It is useful to stress that the solid curve is not a fit 
to these data, but only our best parameterization of the 
SMD in the GOODS-MUSIC sample, which agrees well 
with that of Fontana et al. 2004 in the K20 region and with 
Drory et al. 2005 on the GOODS-South field. Lower panel: 
the mass density computed only on galaxies more massive 
than \og{M^/Mo) — 11. We remark that the value at 
z ~ 0.5 is affected by an undersampling of massive galax- 
ies because of the small size of the volume element. The 
solid line shows the mass density computed with our fidu- 
cial Schechter fit. See text and Tab^f 01, details. The two 
dotted lines represent an exponentially declining evolution 
of the stellar mass density, with timescales t6, 0.6Gyrs ob- 
tained by fitting the data at z < 1.5 and z > 1.5, respec- 
tively. 



4. The comparison with theoretical A-CDM 
models 

In this section we compare our data with a set of recent 
theoretical predictions. In particular, we have included in 
the comparison three recent semianalytic models that in- 



clude the feedback from AGNs on galaxy formation, al- 
beit with different recipes. The first is the latest rendition 
of the semianalytical Durham model (Bow er et al. 2 006 
B06 hereafter), where the feedback from AGNs is ignited 
by the continuous accretion of gas on the central black 
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< z > Schechter Observed Schechter Observed 

8 < log (M/M ) < 13 log (M/Mq) > 11 log (M/M Q ) > 11 

0.50 8.46 8.32+°;^ 8.10 7.89^;^ 

0.70 8.37 8.53±g;g^ 8.01 8.24+.°;o7 

0.90 8.29 8.16to;o3 7.93 7.90±g;^ 

1.15 8.18 8.26to'o2 7.83 7.94l°;°® 

1.45 8.07 7.96t°'o3 7.70 7.57lo'o9 

1.80 7.94 7.90to.o4 7.54 7.50±g:g| 

2.50 7.68 7.60+oo4 7.19 7.16 + .o'o8 

3.50 7.27 7.23±g;^ 6.48 6.60±g;^ 
4.50 7.73' 6.44 
5.50 7.84 7.71 



The Stellar Mass Density at different redshift bins: the second column shows the fitted SMD (integrated from log(M„/Mo) = 8 
to log(M„/Mo) = 13), the third the observed SMD, while the fourth and the fifth columns are the fitted and observed SMD 
for log(M«/M Q ) > 11. 



hole. Such an implementation is conceptually different 
from the M enci et al. 20061 model (M06 in the following), 
where the feedback from AGN is explicitly due to the 
blast wave originated during the active (luminous) phase 
of AGN activity. Differences between these two models 
are better detailed in IMenci et al. 20061 Finally, in the 
semianalytic model of IMonaco et al. 20061 (MFT06) cool- 
ing, infall, star formation, feedback, galactic winds and 
accretion onto black holes are described with a set of 
new simplified models that take into account the multi- 
phase nature of the ISM and its energetics. In particular, 
the quenching of late cooling flows results from the in- 
jection of energy from massive black holes that accrete 
slowly (in Eddington terms) from the cooling flow itself. 
A completely different approach is provided by the hy- 
drodynamical simulations in a cosmological context of 
Nag amine et al. 2005a| |Nagamine et al. 2005b| (and ref- 
erences therein), which have been obtained either with 
a Eulerian mesh code with Total Variation Diminishing 
( |Ryu et al. 2003| ) shock capturing scheme (N-TVD in the 
following) or with a SPH "entropy formulation" method 
dSpringel et al. 2002| N-SPH in the following). Both mod- 
els include radiative cooling and heating, uniform UV 
background, supernova feedback and standard recipes for 
star-formation. In addition, the SPH simulation also in- 
clude the effects of feedback by galactic winds and a mul- 
tiphase ISM, that provides a more accurate modelling of 
the star-formation process. 

We first compare these models with one main result 
of our analysis, the evidence that the evolution of massive 
galaxies is relatively mild up to z ~ 1.5 and then remark- 
ably faster. At this purpose, we plot again in Fig|5j the 
evolution of the stellar mass density in massive galaxies 
as a function of redshift, normalized to the local one, both 
in our data and in the models. 

Since we want to focus on the qualitative behaviour of 
both model and data, and in order to minimize the im- 
pact of different IMFs, we have normalized each model to 



its predicted value at z = 0, and for the models with a 
Kcnnicutt IMF we have used the corresponding charac- 
teristic mass M* at z = 0. 

FigOU clearly shows that all models reproduce the ob- 
served trend in the evolution of the stellar mass density, 
with a slow, steady decay up to z ~ 1.5 and a much faster 
decay thereafter. For all models, the formation of massive 
galaxies is therefore occurring at a fast pace at z > 1.5 — 2, 
and at slower rate at lower redshift, in broad agreement 
with the observed evolution. It has to be remarked, how- 
ever, that the fraction of stellar mass assembled in the 
high redshift phase varies by a factor of two among the 
different models, on which we shall comment later. 

The other important result of our analysis is that 
there is a significant evolution in the density of low mass 
galaxies, that at z ~ 1 are about 40% of the present- 
day number, and that the GSMF is remarkably flat up 
to z ~ 1 — 1.3, with a small evolution of the slope 
(Aa = —0.1) with respect to the local value. A compar- 
ison with the models for this specific feature is provided 
in the upper panel of Fig llOl where we plot the GSMF 
in a broad redshift bin around z ~ 1. It is shown that 
all the models have an unsatisfactory fit to the data, on 
the faint side. Most of them fail to fit the relatively large 
evolution in the density of low mass galaxies, with the re- 
sult of overpredicting their number. In addition, some of 
them, and in particular the M06 one, predict a slope of 
the GSMF much steeper than observed. All these results 
confirm that the observed flatness of the GSMF and most 
of all the significant density evolution at intermediate red- 
shifts is a critical feature, difficult to reproduce by current 
theoretical models. 

We also notice that, always at z ~ 1, most of the mod- 
els roughly predict a correct number of massive galaxies, 
with the exception of the two hydrodynamical simulations. 
The consequence of this comparison is that the ratio be- 
tween massive and low mass galaxies (whose evolution pro- 
vides the "downsizing" scenario) is not reproduced by any 
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Fig. 9. Comparison between the observed Stellar Mass 
Density in massive galaxies (defined as those above M* > 
10 n M©) and the corresponding predictions by recent the- 
oretical models, as a function of redshift. Both data and 
models have been normalized to their own correspond- 
ing local value, for the relevant IMF. The shaded area 
represent the data obtained in the present work, normal- 
ized with two different choices of the local mass func- 
tion (see text for details). The other curves show the 
prediction of theoretical models, as labelled in the leg- 
end: IBower et al. 20061 (B06) , IMenci et al. 20061 (M 06), 
IMonaco et al. 20061 (MFT06), |Nagamine et al. 2005al (N- 
TVD), INagamine et al. 2005 b| (N-SPH). 



model: it is interesting to remark that such a failure is due 
to the overprediction of low mass galaxies, and not to an 
underprediction of massive galaxies. 

In the bottom panels of Fig^| we also compare the 
theoretical GSMFs to the observed ones in the highest red- 
shift bins. The comparison in the redshift range 2 < z < 3 
is probably the most statistically meaningful, albeit any- 
way limited to the higher mass regime, and we briefly 
concentrate on it. It is shown that, at variance with the 
overall agreement found at z ~ 1, the considered model 
span a wide range of predicted densities. Typically, hydro- 
dynamical models tend to overpredict the observed data, 
while the semianalytical renditions appear to be closer 
to the data. In the detail, those that include a "QSO- 
like" feedback from AGNs (M06 and MFT06) tend to un- 
derpredict the data, and those with a "radio-like" AGN 
feedback (B06) tend to overpredict the data. While the 
disagreement is apparently large (by nearly an order of 
magnitude) as far as the number densities are concerned, 
it must be remarked that the effect is emphasized by the 
steep slope of the exponential tail of GSMF, since the 
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Fig. 10. Comparison between the observed GSMF and a 
set of recent theoretical models. As in FigQJ large circles 
represent the mass-complete, Ks-selected sample. Lines 
refer to different models, as in Fig|5] The dashed areas rep- 
resent the confidence regions of the models by Nagamine 
et al 2005ab, where the limits are the predicted mass func- 
tions calculated at z = 2, z — 3 and z = 4. 



offset in mass between the theoretical predictions and the 
data is within a factor of two. Since hydrodynamical mod- 
els do not include feedback from AGNs, while "QSO-like" 
feedback from AGNs is particularly efficient in massive 
halos at high redshift, it is tempting to ascribe the differ- 
ence to this specific feature, although many details of the 
galaxy formation models, especially the cooling and infall 
of gas in the infall-dominatcd halos at high redshift, may 
play a fundamental role as well. 



5. Massive galaxies at z > 4 

As we have described in the previous sections, our sta- 
tistical analysis is limited to z < 4 to take into full ac- 
count the various effects of incompleteness in the sample. 
However, as clearly shown in FigOU our sample includes 
several galaxies at z > 4. Some of these are "drop-out 
galaxies", that are brighter than m z — 26 and detected 
in Ks, the two requirements that we adopted in plotting 
FigCD As such, our sample does not include many other 
"I drop-out galaxies", likely at z > 5, that are detected 
in the GOODS images, since they are either fainter than 
m z — 26 or undetected in Ks. These objects do not make 
up a complete sample and are therefore excluded from our 
current analysis. 

Most notably, however, we have a few very massive 
candidate galaxies at z > 4 in our sample. It is inter- 
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esting and worthwhile to carry on a detailed analysis of 
these objects, especially because they give a significant 
contribution to the mass density of galaxies. In fact, when 
we estimate the stellar mass density for the galaxies with 
z > 4 (FigJSJ) , we note that it is very high, despite it is a 
lower limit: this is more evident in the lower panel of the 
figure where only the mass density of the very massive 
galaxies is shown. Few very massive galaxies are respon- 
sible for such high value of the mass density: at z = 3 — 4 
two massive galaxies give 34% of the total mass density, 
while at z = 5 — 6 the mass density of only three objects 
is 94% of the total value. 

It is interesting to investigate what kind of objects are 
these. If we select the galaxies with log (M„,/M©) > 11, 
we find three objects, as shown in Fig. [3J two of them are 
Ks and Zgso-selected, and one is only Ks-selected. These 
objects are typically quite red for their redshifts, with a 
rest-frame U — V ~ 1.5. The best-fitting templates for 
these objects are either passively evolving models, char- 
acterized by a very short timescale for star formation t 
with E(B — V) ~ or by a constant star formation model 
with large amount of dust, E(B — V) > 0.5. Two of these 
objects are displayed in the left panel of Fig. ^2 It is im- 
portant to stress that the redshift determinations of all 
these objects are very uncertain, as it is shown by the % 2 
as a function of redshift in the same figure. In practice, on 
the basis of the available data one can only conclude that 
these objects are most likely at z > 2, but a sound red- 
shift determination is not available. In particular, other 
redshift solutions can be find with a similar probability at 
z = 2-3. 

We note that these objects might resemble the object 
presented in Mobasher et al. 2005, that is also a candidate 
massive, passively evolving galaxy at z ~ 6.5. We have 
also extracted from our raw catalog the photometry of 
the Mobasher et al. (2005) object (that is slightly fainter 
than our limiting magnitude in Ks), to check how it is 
classified with our photometry, that is based on a different 
image set (Mobasher et al. (2005) used a combination of 
UDF images taken with ACS and NICMOS). With our 
photometry, such an object is quite similar to our z ~ 2 
candidates in Fig. 1111 with a flat probability distribution 
for the photometric redshift from z = 2 to 7 and the best 
fit solution at z = 2.25 with Eb-v = 1-1: in this case the 
stellar mass turns out to be 2 ■ 10 11 M Q , quite typical at 
this redshift. 

At this point, it is interesting to seek for similar galax- 
ies in the GOODS-MUSIC sample, i.e. galaxies that have 
a comparable spectral distribution, best-fitting template 
around z = 2 — 3 and with a probability distribution ex- 
tending up to z ~ 5 — 6. Such objects indeed exist in our 
sample: we find in total three objects, a number compara- 
ble to the number of the massive object at higher redshift. 
Two of them are shown in the right panel of Fig. ^2 The 
fraction of objects with uncertain redshifts is negligible at 
z = 2 — 3 where we have estimated MF and mass density, 
such that we do not expect that they introduce a signi- 
ficative uncertainty in the statistical analysis. In any case, 




MA) MA) 



Fig. 11. Left: Spectral Energy Distributions of two galax- 
ies with photometric redshift > 4, that would result very 
massive M* ~ 5 x 10 1 M©. In the inset, the x 2 as a func- 
tion of redshift of the photometric redshift is shown. Right 
panel: Spectral Energy Distributions of two objects with 
similar colours and uncertainty in redshift, that are as- 
signed a photometric redshift ~ 2. 



the effect of such uncertainty is taken into account by the 
Monte Carlo estimate of the error budget in the GSMF. 

We can draw two different conclusions from this exer- 
cise. First, although we cannot exclude that some of these 
objects are actually very massive galaxies at z > 4, we 
have to await for a more robust determination of their 
redshift before including them in a firm estimate of the 
stellar mass density at high redshift. 

Second, the properties of the GSMF can give some 
hints on the nature of these objects. Probably the most 
important result of our analysis is that we can exclude 
that these objects are typical in a statistical sense. Our 
sample is indeed complete for passively evolving galaxies 
of stellar mass M* > 10 n M© up to z ~ 4, and we can 
firmly estimate a decrease of at least a factor of ten in 
their density at such high redshift. A scenario where most 
of present-day early galaxies formed at even higher red- 
shift z > 4 and evolved passively thereafter is therefore 
excluded by our analysis, such that we can conclude that 
the large mass density contained in passively evolving sys- 
tems that would arise from our z > 4 candidates is either a 
statistical fluctuation, or results from incorrectly assigned 
photometric redshifts. Future analysis will hopefully clar- 
ify this point. 
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6. Summary and Discussion 

In this work we have presented an analysis of the evo- 
lution of the Galaxy Stellar Mass Function (GSMF) and 
of the corresponding stellar mass density up to z ~ 4. It 
has been obtained from the GOODS-MUSIC sample, a 
-fTs-selected catalog of 2931 galaxies with 14 bands pho- 
tometry, extending from 0.35 to 8/im, extracted from the 
public data of the GOODS-South survey. We have derived 
accurate stellar masses for this large sample of galaxies, 
adopting the standard technique of fitting spectral synthe- 
sis models (Bruzual & Chariot 2003) on the multiwave- 
length data. To compare with previous works, we have 
used the Salpeter IMF, ranging over a set of metallici- 
ties (from Z = 0.02Z Q to Z = 2.5-Z© and dust extinc- 
tion (0 < E(B — V) < 1.1, with a Calzetti extinction 
curve). On this data set, we have computed the Galaxy 
Stellar Mass Function and the resulting stellar mass den- 
sity with the usual 1/V max formalism. To provide a vi- 
sualization of the GSMF free of the well known biases 
of the l/V max approach (namely, the sensitivity to large 
scale structures and the arbitrarity of the binning) , we also 
provide a global fit to the data with a Schechter function 
with smooth, redshift -dependent free parameters ("STY 
approach"). 

We have carefully discussed the selection effects that 
play a role in our K s-selected sample, and that must be 
taken into account in the analysis. We have shown that, 
while we can detect massive (M» > 10 11 M Q ) passively 
evolving galaxies up to z ~ 4, our sample becomes pro- 
gressively biased against star-forming, dust-enshrouded 
objects with E(B — V) > 0.5 already at z > 2. Since 
there are evidences that these objects do exist at z < 2, it 
is possible that our census of galaxies at z > 2 is incom- 
plete, and that our estimates in this redshift range (as in 
any if s-selected sample) should be considered as a lower 
limit. 

The major results that we have found with these data 
are the following: 

- We show that the inclusion of the Spitzer data (from 
3.5 to 8 fjxa) significantly improves the reliability of the 
mass estimate at z > 2, as expected. At lower redshift, 
there is no significant improvement and the observed scat- 
ter is consistent with that induced by the model degener- 
acy. 

- We confirm the well known trend of global decline 
of the stellar mass density with redshift. The total mass 
density is about a factor 2x lower than local at z ~ 1, and 
about 10% of the local at z ~ 3. 

- We compute the GSMF in several redshift bins, from 
z = OA to z — A, and we show that a simple scaling of 
the Schechter parameters is able to provide a smoothly 
evolving rendition of the GSMF. 

- Restricting only to galaxies with M« > 10 11 M Q , their 
mass density evolves relatively mildly up to z ~ 1 — 1.5. At 
z ~ 1, their integrated mass density is about 50% of the 
local value. At z > 1.5, massive galaxies become to evolve 
much faster, such that at z ~ 3 they provide at most one 



tenth of the local density. This trend may be described by 
two exponential e-folding times, i.e. p(z) oc e _t / T , with 
r ~ 6Gyrs at z < 1.5 and r ~ 0.6Gyrs at z > 1.5, although 
this number may be biased low because of the selection 
effects described above. 

- As far as low mass galaxies are concerned, we show 
that the GSMF remains remarkably flat (it steepens by 
less than 0.1 for each unit redshift) up to z ~ 1 — 1.3, the 
highest redshift bin where it can be reliably measured. At 
the same time, a sensible decrease occurs in the number 
of low mass galaxies: the density of galaxies with M ~ 
10 10 M© at z ~ 1 is 4 times lower than the local one. 

- In this context, we finally show that there is a clear 
evidence for a differential evolution of the Galaxy Stellar 
Mass Function up to z ~ 1.5, with less massive galaxies 
evolving more than massive ones. Such a trend is evident 
in FigO] that shows that the GSMF in the redshift bins 
up to z — 1.3 remains progressively closer to the local 
one for increasing stellar masses, and it is also substan- 
tiated by two quantitative evidences. First, from z = 
to z ~ 1.5 the total stellar mass density decreases more 
than that due to massive galaxies only (i.e. those with 
M* > 10 U M Q ). More intriguingly, the STY fit to the 
overall evolution of the GSMF shows that the characteris- 
tic mass M* increases up to z ~ 1.3, reflecting the change 
in the shape of the GSMF, that from z ~ to z ~ 1.3 is 
progressively more skewed toward higher masses. 

- After z ~ 1.5, however, the evolution of the mass 
density contained in massive galaxies is much faster, as 
shown by Fig0 showing that this is their main epoch of 
formation. 

It is straightforward to see that such differen- 
tial evolution of the GSMF is a natural conse- 
quence of the "downsizing" scenario for galaxy evo- 
lution, that has been found in many independent 
surveys. Following the original definition by Cowie 
et al. 1996, it has been shown by several au- 
thors fe.g. IBrinchmann" fc Ellis 20001 IFontana et al. 20031 
IPerez-Gonzales et al 20051 IFeulner et al. 2005)) that the 
specific star-formation rate increases with redshift, show- 
ing that massive galaxies become progressively more ac- 
tively star-forming as redshift increases. By looking at the 
GSMF, we observe the consequences of such a trend: more 
massive systems form in a vigorous phase at high redshift, 
that is largely complete at z ~ 1 — 1.5, such that the corre- 
sponding section of the GSMF is already close to the local 
value. Lower mass systems, on the contrary, continue to 
grow their stellar content at even lower redshifts, such that 
the increase of the faint side of the GSMF is large also at 
low z. 

We remark that this picture is not plagued by the bi- 
ases against the detection of dusty massive galaxies at high 
z that exist in our sample, as in any Ks complete one. As 
we discussed in the text, we might be missing high red- 
shift, star-forming dusty galaxies with large extinction, 
but we are essentially complete with respect to passively 
evolving, dust-free galaxies of M* > 10 11 M Q up to z ~ 4. 
The low density of high redshift massive systems provides 
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decisive evidence that the formation and assembly of lo- 
cal, massive bulges and ellipticals did not form in a single 
phase at very high redshift. 

We have explored whether such a picture is qualita- 
tively consistent with the predictions of most recent the- 
oretical models for galaxy formation. We have compared 
our data with a set of theoretical models, including semi- 
analytic models (JBower et al. 20061 IMenci et al. 2 006 
IMonaco et al. 2006jl as well as hydrodynamical ones 
( |Nagamine et al. 200"5a| |Nagamine et al. 2005b| ). All 
these models predict a relatively mild evolution of the 
stellar mass density contained in massive galaxies from 
z = to z — 1.5, that is broadly consistent with the 
observed data, and a much faster evolution at higher z, 
again in agreement with the data. 

For all the models, however, at intermediate redshifts 
the match with the slope and the normalization of the 
GSMF at intermediate or low masses (M* ~ 10 10 Af©) is 
still critical, when not even poor. Interestingly, this implies 
that the "downsizing" scenario (that is based on the evo- 
lution ratio between massive and low mass galaxies) is not 
reproduced by any model because of the overprediction of 
low mass galaxies, and not because of an underprediction 
of the massive ones. 

At high redshift, the detailed agreement with the ob- 
served data is still sensitive to the description of the phys- 
ical processes inserted in the models. Hydrodynamical 
models tend to overpredict the observed mass density, 
as already noted by the authors (Nagaminc et al. 2004| ), 
while the semianalytic models that include feedback from 
AGNs are closer to the data. 

In conclusion, the new data presented in this work pro- 
vide an overall description of the rise of the stellar mass, 
and in particular of that residing in massive galaxies, in 
which about half of such stellar mass appear to have been 
assembled during the first 2-4 Gyrs after recombination, 
followed by a milder increase over the remaining cosmic 
time. Although encouraging, the comparison with the the- 
oretical expectations provides evidences that some funda- 
mental physical processes, likely affecting both low and 
high mass galaxies, are still incorrectly represented, and 
that the density of massive galaxies at high redshift is 
indeed a very useful test for these models. 

At this purpose, more effort is needed to improve the 
reliability of the observational estimates at intermediate 
and high redshift. From the observational point of view, 
larger spectrophotometric surveys on independent fields 
are definitely necessary to smear out the effects of sample 
variance. These might in particular affect two of our major 
results, namely the increase of the characteristic mass M* 
up to z ~ 1 and the amount of massive galaxies at z > 2. 
In addition to such (obvious) caveats, several systematic 
effects are still to be properly minimized. To mention a 
few, those related to the choice of the IMF and to the 
differences that may arise from the treatment of post- 
AGB stars in spectral synthesis models (|Maraston 2005|l . 
and the possible contribution of dust-enshrouded galaxies 
to the overall mass budget. 



In addition, these findings rise important questions 
about the physical processes that led to the rise of the 
stellar mass density in massive galaxies: what is the phys- 
ical nature of the galaxies that contribute to the stellar 
mass density at high redshifts, and what is the physi- 
cal mechanism that drove this rise, i.e. how much it is 
related to star-formation occurring within the observed 
galaxies as compared to the contribution from merging 
processes. Although present-day surveys are starting to 
explore these issues, these will remain among the more 
challenging questions of the next years. 
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